home *** CD-ROM | disk | FTP | other *** search
/ NeXT Education Software Sampler 1992 Fall / NeXT Education Software Sampler 1992 Fall.iso / Programming / Source / 2DLab / geometry.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-01-14  |  3.9 KB  |  190 lines
  1. #include <math.h>
  2. #include "voronoi.h"
  3.  
  4. geominit()
  5. {
  6.     float sn;
  7.  
  8.     freeinit (&efl, sizeof (struct Edge));
  9.     nvertices = 0;
  10.     nedges = 0;
  11.     sn = nsites + 4;
  12.     sqrt_nsites = sqrt(sn);
  13.     deltay = ymax - ymin;
  14.     deltax = xmax - xmin;
  15.     siteidx = 0;
  16. }
  17.  
  18. struct Edge *bisect(s1,s2)
  19. struct    Site *s1,*s2;
  20. {
  21.     float dx,dy,adx,ady;
  22.     struct Edge *newedge;
  23.  
  24.     newedge = (struct Edge *) getfree(&efl);
  25.  
  26.     newedge->reg[0] = s1;
  27.     newedge->reg[1] = s2;
  28.     ref(s1); 
  29.     ref(s2);
  30.     newedge->ep[0] = (struct Site *) NULL;
  31.     newedge->ep[1] = (struct Site *) NULL;
  32.  
  33.     dx = s2->coord.x - s1->coord.x;
  34.     dy = s2->coord.y - s1->coord.y;
  35.     adx = dx>0 ? dx : -dx;
  36.     ady = dy>0 ? dy : -dy;
  37.     newedge->c = s1->coord.x * dx + s1->coord.y * dy + (dx*dx + dy*dy)*0.5;
  38.     if (adx>ady) {    
  39.     newedge->a = 1.0; 
  40.     newedge->b = dy/dx; 
  41.     newedge->c /= dx;
  42.     } else {    
  43.     newedge->b = 1.0; 
  44.     newedge->a = dx/dy; 
  45.     newedge->c /= dy;
  46.     }
  47.     newedge->edgenbr = nedges;
  48.     out_bisector(newedge);
  49.     nedges += 1;
  50.     return(newedge);
  51. }
  52.  
  53. struct Site *intersect(el1, el2, p)
  54. struct Halfedge *el1, *el2;
  55. struct Point *p;
  56. {
  57.     struct Edge *e1,*e2, *e;
  58.     struct Halfedge *el;
  59.     float d, xint, yint;
  60.     int right_of_site;
  61.     struct Site *v;
  62.  
  63.     e1 = el1->ELedge;
  64.     e2 = el2->ELedge;
  65.     if (e1 == (struct Edge*)NULL || e2 == (struct Edge*)NULL) 
  66.     return ((struct Site *) NULL);
  67.     if (e1->reg[1] == e2->reg[1]) 
  68.     return ((struct Site *) NULL);
  69.     d = e1->a * e2->b - e1->b * e2->a;
  70.     if (-1.0e-10<d && d<1.0e-10) 
  71.     return ((struct Site *) NULL);
  72.     xint = (e1->c*e2->b - e2->c*e1->b)/d;
  73.     yint = (e2->c*e1->a - e1->c*e2->a)/d;
  74.     if ((e1->reg[1]->coord.y < e2->reg[1]->coord.y) ||
  75.         (e1->reg[1]->coord.y == e2->reg[1]->coord.y &&
  76.         e1->reg[1]->coord.x < e2->reg[1]->coord.x) ) {    
  77.     el = el1; 
  78.     e = e1;
  79.     } else {    
  80.     el = el2; 
  81.     e = e2;
  82.     }
  83.     right_of_site = xint >= e->reg[1]->coord.x;
  84.     if ((right_of_site && el->ELpm == le) ||
  85.        (!right_of_site && el->ELpm == re)) 
  86.     return ((struct Site *) NULL);
  87.     v = (struct Site *) getfree(&sfl);
  88.     v->refcnt = 0;
  89.     v->coord.x = xint;
  90.     v->coord.y = yint;
  91.     return(v);
  92. }
  93.  
  94. /* returns 1 if p is to right of halfedge e */
  95. int right_of(el, p)
  96. struct Halfedge *el;
  97. struct Point *p;
  98. {
  99.     struct Edge *e;
  100.     struct Site *topsite;
  101.     int right_of_site, above, fast;
  102.     float dxp, dyp, dxs, t1, t2, t3, yl;
  103.  
  104.     e = el->ELedge;
  105.     topsite = e->reg[1];
  106.     right_of_site = p->x > topsite->coord.x;
  107.     if (right_of_site && el->ELpm == le) 
  108.     return(1);
  109.     if (!right_of_site && el->ELpm == re) 
  110.     return (0);
  111.     if (e->a == 1.0) {    
  112.     dyp = p->y - topsite->coord.y;
  113.     dxp = p->x - topsite->coord.x;
  114.     fast = 0;
  115.     if ((!right_of_site && e->b<0.0) || (right_of_site && e->b>=0.0) ) {    
  116.         above = dyp>= e->b*dxp;    
  117.         fast = above;
  118.     } else {    
  119.         above = p->x + p->y*e->b > e->c;
  120.         if (e->b<0.0) 
  121.         above = !above;
  122.         if (!above) 
  123.         fast = 1;
  124.     }
  125.     if (!fast) {    
  126.         dxs = topsite->coord.x - (e->reg[0])->coord.x;
  127.         above = e->b * (dxp*dxp - dyp*dyp) <
  128.                     dxs*dyp*(1.0+2.0*dxp/dxs + e->b*e->b);
  129.         if (e->b<0.0) 
  130.         above = !above;
  131.     }
  132.     } else {    /*e->b==1.0 */
  133.     yl = e->c - e->a*p->x;
  134.     t1 = p->y - yl;
  135.     t2 = p->x - topsite->coord.x;
  136.     t3 = yl - topsite->coord.y;
  137.     above = t1*t1 > t2*t2 + t3*t3;
  138.     }
  139.     return (el->ELpm==le ? above : !above);
  140. }
  141.  
  142. void
  143. v_endpoint(e, lr, s)
  144. struct Edge *e;
  145. int    lr;
  146. struct Site *s;
  147. {
  148.     e->ep[lr] = s;
  149.     ref(s);
  150.     if (e->ep[re-lr]==(struct Site *)NULL) 
  151.     return;
  152.     out_ep(e);
  153.     deref(e->reg[le]);
  154.     deref(e->reg[re]);
  155.     makefree(e, &efl);
  156. }
  157.  
  158. float dist(s,t)
  159. struct Site *s,*t;
  160. {
  161.     float dx,dy;
  162.     dx = s->coord.x - t->coord.x;
  163.     dy = s->coord.y - t->coord.y;
  164.     return(sqrt(dx*dx + dy*dy));
  165. }
  166.  
  167.  
  168. int makevertex(v)
  169. struct Site *v;
  170. {
  171.     v->sitenbr = nvertices;
  172.     nvertices += 1;
  173.     out_vertex(v);
  174. }
  175.  
  176.  
  177. deref(v)
  178. struct    Site *v;
  179. {
  180.     v->refcnt -= 1;
  181.     if (v->refcnt == 0) 
  182.     makefree(v, &sfl);
  183. }
  184.  
  185. ref(v)
  186. struct Site *v;
  187. {
  188.     v->refcnt += 1;
  189. }
  190.